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TSP - TREND SURFACE PLOTTING 


DECUS Program Library Write-up 


DECUS No. 8-324 


ABSTRACT 


This program consists of a group of segments for which the output of one segment is the input 
of the next. It enables the significance of linear, quadratic, and cubic trend surfaces to 
be determined for each of a number of variables on the co-ordinates of their distribution in 
two-dimensional space. The significance of the various surfaces having been determined, 
the fitted surface can be plotted by means of a contour map on the teletype. 


TAPES REQUIRED 


1. Form of program tapes - All program tapes are written in the PDP-8 FORTRAN-D 
language and are in the source language for ease of modification. There are seven segments 
in the complete trend surface plotting program, as follows: 


(a) TSP2 - program to compute the regressor matrices from the original data tape of the grid 
co-ordinates; 


(b) TSP3 - program to invert the regressor matrices computed by TSP2; 

(c) TSP4 - program to compute sums of squares and products vectors for dependent variables; 

(d) TSP5 - program to compute regressions of dependent variables on the linear, quadratic, 
and cubic expressions of the grid co-ordinates; 

(e) TSP6 - program to compute the significance of the computed regressions; 


(f) TSP7 - program to compute trend surface parameters. (N.B. this program will usually 
require to be modified to suit individual applications. For advice on how these modifications 
should be made, consult the author.) 


(g) TSP8 - program to plot trend surfaces. 

2. Form of data tapes - Three data tapes are required for this program: 

(a) Data tape fot TSP2 - This tape should consist of the number of sample points, followed by 
the grid co-ordinates for each point, e.g. 


41 

3.44 7.29 

3.30 7.50 

3.15 7.70 

3.30 7.70 etc. 








(b) Data tape for TSP4 - This tape should consist of the successive values of the dependent 
variable, with the sampling points in the same order as the data tape for TSP2. As many 
dependent variables as are required can be included on this data tape, eash as a single column 
vector. 

(c) Data tape for TSP5 - A short data tape containing the means of the linear, quadratic, 
and cubic terms of the grid co-ordinates is required, and can be derived from the printed 
results of TSP2. 

OPERATING INSTRUCTIONS 

All of the program segments follow the usual operating procedures for the PDP-8 disk operating 
system. The sequence of the various calculations is given in Figure 1 . Detailed instructions 
for the individual segments are as follows: 

(a) TSP2 - The data tape for this segment is placed in the high-speed tape reader before 
continuing after the teletype has printed "READY. " The high-speed punch should be switched 
on while a summary of the data is being typed. The output from this segment should be saved, 
and used as the input to TSP3. 

(b) TSP3 - The output from TSP2 should be placed in the high-speed reader, and the high¬ 
speed punch switched on before continuing after "READY." The output from this segment 
should be saved and used as an input for TSP5. 

(c) TSP4 - The data tape for TSP4 (see 2 (b) above) should be placed in the low-speed 
tape reader, and the data tape for TSP2 in the high-speed reader, and the high-speed 
punch switched on before continuing after "READY. " The output from this segment should 
be saved and used as an input to TSP5. 

(d) TSP5 - The data tape containing the means of the linear, quadratic, and cubic terms of 
the grid co-ordinates (see 2 (c) above) should be placed in the high-speed reader before 
continuing after "READY." The program will then pause for the insertion of the output from 
TSP4 in the slow-speed reader and the output from TSP3 in the high-speed reader, before 
continuing. After computing the regressions for the first dependent variable, the program will 
pause for the output from TSP3 to be replaced in the high-speed reader before continuing to 
read the next dependent variable. 

(e) TSP6 - Before running the TSP6 segment, it is necessary to create a data tape from the 
printed results of TSP5. This is done by punching the number of sample points, followed by 
the total sum of squares, sum of squares due to regression, and coefficient of determination 
for the linear, quadratic and cubic regressions. This data tape should be placed in the high¬ 
speed reader before continuing after "READY. " 

(f) TSP7 - On continuing after "READY, " the program will pause for the entry of the re¬ 
gression coefficients of the significant trend surface. Note that the values for all nine co¬ 
efficients must be entered even if they are zero. The coefficients must then be followed by 
the appropriate constant term, and by the starting value for the contours to be drawn and the 
contour interval. As a rough guide, the starting value should be a little below the minimum 
value recorded in the dependent variable and the contour interval about one third of the 
standard deviation of the dependent variable. The output from TSP7 is stored directly on the 
disk, in preparation for the use of TSP8. 
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(g) TSP8 - Before continuing after "READY, " it is advisable to title the trend surface plot 
with the teletype switched to "local, " followed by five or six line feeds. If an output tape 
is required, switch on the low-speed printer before continuing. 

OUTPUT 

The output from the various program segments is as follows: 

(a) TSP2 - The means and standard deviations of the linear, quadratic, and cubic terms are 
typed, and the corrected sums of squares and products for the linear, quadratic, and cubic 
regressions are output on the high-speed punch in a form suitable for direct re-input to TSP3. 

(b) TSP3 - The inverse regressor matrices are output on the high-speed punch in a form 
suitable for direct re-input to TSP5. 

(c) TSP4 - This program outputs a vector of sums of squares and products for each dependent 
variable, in a form suitable for direct input to TSP5. 

(d) TSP5 - This program prints the number of sampling points, the partial regression coefficients 
and constant, the total sum of squares and the sum of squares due to regression and the pro¬ 
portion of the total variability accounted for by the regression for the linear, quadratic, and 
cubic trend surfaces. 

(e) TSP6 - The program prints the degrees of freedom and variance ratio testing the signifi¬ 
cance of the linear, quadratic, and cubic components of the trend surface. 

(f) TSP7 - The parameters for the printing of a contour map of the trend surface are stored 
directly on the disk. 

(g) TSP8 - The contour map of the fitted trend surface is typed. An example of a typical 
map is given in Figure 2. 

STORAGE AND LIMITATIONS 

Normal for FORTRAN-D 

There are no limits on the number of points that can be used in the computation of the trend 
surfaces . The scale of the map that can be printed will usually be limited by the width of 
the printed line on the teletype. 

METHOD 

The method follows closely that described by the following papers: 

J. W. Harbaugh, BALGOL program for trend-surface mapping. Special Distribution Publica¬ 
tion, University of Kansas, 1963. 

M. O'Leary, R. H. Lippert, and O. T. Spitz, FORTRAN IV and map program for computation 
and plotting of trend surfaces for degrees 1 through 6, Computer Contribution 3, State Geologi¬ 
cal Survey, University of Kansas, 1966. 



Figure 1 


Flow diagram for segments of trend-surface plotting program 
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Figure 2 . 
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♦OFT-S 
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PROGRAM TO COMPUTE REGRESSOR MATRICES TSP2 
DIMENSION SX(9), 55(81), X(9) 

READ 2,101, N 
FORMAT (I) 

DO 10 1=1,9 

SX(l)=0.0 

CONTINUE 

DO 70 J=1,81 

SS(J)=0.0 

CONTINUE 

DO 16 1=1, N 

READ 2,105,X(1),X(2) 

FORMAT (E, E) 

X(3)=X(1)*X(1) 

X(4)=X(2)*X(2) 

X(5)=X(1)*X(2) 

X(6)=X(3)*X(1) 

X(7)=X(4)*X(2) 

X(8)=X(3)*X(2) 

X(9)=X(4)*X(1) 

E 1=1 

DO 20 J=1,9 
DO 20 K=1,9 
KK=J+9*(K-1) 

SS(KK)=SS(KK)+(X(J)-SX(J) )*(X(K)-SX(K) )*(1 . 0 -] .0/EI) 

CONTINUE 

DO 16 J=1,9 

SX(J)=SX(J)+(X(J)-SX(J) )/EI 

CONTINUE 

EN=N 

DO 35 1=1,9 
K=l+9*(l-l) 

X(I)=SQTF(SS(K)/(EN-1 . 0 )) 

CONTINUE 
TYPE 101, N 
DO 103 1=1,9 
TYPE 104, SX(I), X(l) 

CONTINUE 
FORMAT (/, E, E) 

M=2 

WRITE 2,101 ,M 

WRITE 2,105, SS(1), SS(2), SS(10), SS(11) 

M=5 

WRITE 2,101 ,M 
DO 40 1=1,37,9 
DO 40 J=1,5 
K=H-J-1 

WRITE 2,105, SS(K) 



CONTINUE 

M=9 

WRITE 2, 1J01 ,M 

DO 60 1=1,81 

WRITE 2,1J05,SS(I) 

CONTINUE 

STOP 

END 



L 

AM TO INVERT MATRICES TSP3 


230 

DIMENSION A(156) 

READ 2, 3, N 

3 

FORMAT(I) 

5 

LAST =N*N 

11=1 

12=LAST-N+1 

DO 100 J=l, N 

DO 101 1=11,12, N 

READ 2,5, A(l) 

FORMAT (E) 

101 

CONTINUE 

100 

11=11+1 

12=12+1 

CONTINUE 

6 

WRITE 2,6, N 

FORMAT (/, 1) 

149 

DO 150 J=1, N 

105 

DO 105 1=1, N 

A (LAST+l)=0.0 
CONTINUE 

106 

A(LAST+J)=1.0 

PVT=A(J) 

J3=LAST+J 

DO 106 KP=J,J3,N 

A(KP)=A(KP)/PVT 

CONTINUE 

107 

DO 110 KRT=1, N 

IF (KRT-J) 107,110,107 
KR1=KRT 

109 

KR2=KR1+LAST 

KPR=J 

RWC=A(KR1) 

DO 109 KR=KR1, KR2, N 

A(KR)=A(KR)-RWC*A(KPR) 

KPR=KPR+N 

CONTINUE 

110 

CONTINUE 

111 

DO 111 1=1, LAST 
A(I)=A(I+N) 

CONTINUE 

150 

CONTINUE 

200 

ILNE=4 


NR=NRC=1 

11=1 

I2=LAST-N+1 

KONT =1 

DO 210 JP=1, N 
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212 

216 

217 

215 

220 

221 

214 

C 213 

21J0 

END 


DO 213 1=11,12, N 
WRITE 2,212, A(l) 
FORMAT(E) 

IF(N-NRC)216,216,215 

IF(N- NR)213,213,217 

NR=NR+1 

NRC=KONT=l 

ILNE=4 

GO TO 213 

IF(KQNT-ILNE)214,220,214 
ILNE=ILNE+4 
WRITE 2,221 
FORMAT (/,/) 
KONT=KONT+l 
NRC=NRC+1 
CONTINUE 
11 = 11+1 
12 = 12+1 
CONTINUE 
GO TO 23J0 
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C PROGRAM TO COMPUTE DEPENDENT VECTOR TSP 4 

DIMENSION X(l^),SS(ljef), SX(10) 

1 READ 2,101, N 

101 FORMAT (/, I) 

WRITE 2, 101, N 
DO 10 1=1, 10 
X(l) = 0.0 

SS(I) = 0.0 
SX(I) = 0.0 
10 CONTINUE 

DO 20 1=1, N 
READ 2,102, X(I),X(2) 

102 FORMAT (E,E) 

X(3)=X(1)*X(1) 

X(4)=X(2)*X(2) 

X(5)=X(1)*X(2) 

X(6)=X(3)*X(1) 

X(7)=X(4)*X(2) 

X(8)=X(3)*X(2) 

X(9)=X(4)*X(1) 

E 1=1 

READ 1,102,X(10) 

DO 40 J=1,10 

SS(J)=SS(J)+(X(J)-SX(J) )*(X(10)-SX(10) )*(1.0-1.0/EI) 
SX(J)=SX(J)+(X(J)-SX(J) )/EI 
40 CONTINUE 

20 CONTINUE 

EN=N 

DO 50 1=1,10 
WRITE 2, 102, SS(I) 

50 CONTINUE 

WRITE 2, 102, SX(10) 

PAUSE 
GO TO 1 
END 
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C 


60 

4 

] 0 ] 


102 

]0 

3 



20 


04 

30 



50 


103 

5 

6 


PROGRAM TO COMPUTE REGRESSION STATISTICS 
DIMENSION SXY(11), SX(81), B(10) f X(10) 

DO 60 1=1,9 
READ 2, 102,X(\) 

CONTINUE 

PAUSE 

READ 1 , 1J01 , N 
TYPE 101 , N 
FORMAT (/,/, I) 

DO 10 1=1,11 
READ 1,102,SXY(I) 

FORMAT (/, E) 

CONTINUE 
READ 2,1J01 ,M 

DO 20 1=1,MM 
READ 2, 1J02,SX(I) 

CONTINUE 

rss=jM 

K=1 

DO 30 1=1, M 
B(IHM 
DO 40 J=1,M 
B(I)=B(I)+SX(K)*SXY(J) 

K=K+1 

CONTINUE 

CONTINUE 

B(1 ( 0)=SXY(11) 

DO 50 1=1 ,M 
B(lje()=B(l^)-B(l)*X(l) 

TYPE ljef2,B(l) 

RSS=RSS+B(I)*SXY(I) 

CONTINUE 
TYPE 102,B(10) 

RR=RSS/SXY(1^) 

TYPE 103,SXY(lJ0),RSS,RR 
FORMAT (/, E, E, E) 

IF (M-9) 5,6,5 
GO TO 3 
PAUSE 
GO TO 4 
END 
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C PROGRAM TO COMPUTE SIGNIFICANCE OF REGRESSIONS TSP 6 

DIMENSION SS(9) 

READ 2, ]fl ], N 

101 FORMAT (I) 

1 DO 10 1=1,9 

READ 2,102,SS(I) 

102 FORMAT (E) 

10 CONTINUE 

NR=2 

NDF=N-3 

EN=NDF 

RSD=(SS(1)-SS(2) )/EN 
FR=(SS(2)/2.0)/RSD 
TYPE 103, NR, NDF, FR 

103 FORMAT </,/,/,/,/, I, I, E) 

NR=3 

NDF-N-6 
EN=NDF 

RSD=(SS(4)-SS(5) )/EN 
FR=((SS(5)-SS(2) )/3.0)/RSD 
TYPE 104, NR, NDF, FR 

104 FORMAT (/, I, I, E) 

NR=4 

NDF=N-10 

EN=NDF 

RSD=(SS(7)-SS(8) )/EN 
FR=( (SS(8)-SS(5) )/4.0)/RSD 
TYPE 104, NR, NDF, FR 
GOTO 1 
END 


* 
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PROGRAM TO COMPUTE TREND SURFACE PARAMETERS TSP7 
DEFINE DISK 

10 

101 

102 

5 

DIMENSION K(74),B(10) / V(9) 

DO 10 l=U74 

READ 2,101, K(l) 

CONTINUE 

FORMAT (1) 

FORMAT (E) 

N5=l 

20 

N1=0 j 

V24=8.0 

DO 20 1=1,9 

ACCEPT 1102,6(1) 

CONTINUE 

ACCEPT 102, B(10) 

ACCEPT 103, V21,V22 

C 103 

4 

3 

31 

FORMAT (E, E) 

V23=2.05 

IF (K(N5) ) 31,1,1 

N0=-K(N5) 

DO 30 1=1, N|0 

N1=N1+1 

30 

11 

V23=V23+0.05 

WRITE 3,101, K(N5) 

CONTINUE 

N5=N5+1| 

32 

IF (Nl-60) 3,32,3 

V24=V24-0.1 

N1=0 

IF (N5-74) 4,4,33 

33 

STOP 

GO TO 5 

• 

N3=K(N5) 

N4=0 

9 

V(1)=V23 

V(2)=V24 

V(3)=V(1)*V(1) 

V(4)=V(2)*V(2) 

V(5)=V(1)*V(2) 

V(6)=V(3)*V(1) 

V(7)=V(4)*V(2) 

V(8)=V(3)*V(2) 

V(9)=V(4)*V(1) 

N0=0.0 

DO 40 1=1,9 

V0=V(l)*B(l)+V0 
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CONTINUE 

V0=V0+B(10) 

N2=0 

V25=V21 

V25=V25fV22 

N2=N2+1 

IF (V25-V0) 34,34,35 
IF (N2-19) 7,35,35 
WRITE 3,101, N2 
N4=N4+1 
N1=N1+1 
V23=V23+0.05 
IF (N4-N3) 9,11,9 
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32 


20 

99 


12 

m 



1 

101 

2 

102 


3 

03 

4 


5 

104 


6 



7 

105 

8 


9 

106 

10 

11 

107 

21 


13 

108 


TREND SURFACE PLOTTING PROGRAM TSP 8 
DEFINE DISK 
DIMENSION K(60) 

DO 40 1=1,20 
DO 20 J=1,60 
READ 3,99, K(J) 

CONTINUE 

FORMAT(I) 

DO 30 J=1,60 
IF(K(J) ) 21,22,22 
TYPE 100 
FORMAT 
GO TO 30 
TYPE 101 
FORMAT ("0") 

GO TO 30 
TYPE 102 
FORMAT (" ") 

GO TO 30 
TYPE 103 
FORMAT ("1") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 104 
FORMAT ("2") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 105 
FORMAT ("3") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 106 
FORMAT ("4") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 107 
FORMAT ("5") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 108 
FORMAT ("6") 

GO TO 30 
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14 

15 
109 

61 

17 

11 

81 

91 

111 

22 

03 

121 

04 

131 


TYPE 102 
GO TO 30 
TYPE 109 
FORMAT ("7") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 110 
FORMAT ("8") 

GO TO 30 
TYPE 102 
GO TO 30 
TYPE 111 
FORMAT ("9") 

GO TO 30 
LIST=K(J) 

GO TO (1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19),LIST • 

CONTINUE 

TYPE 112 

FORMATS) 

CONTINUE 
TYPE 113 

FORMAT (/,/,/,/,/,/,/,/,/,/ ) 

PAUSE 
GO TO 23 
END 


* 
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